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Abstract 


The incompressible and steady flow of fluid over a square cylinder 
represents an idealization of many industrially important processes, e.g. 
chimneys, tubular heat exchangers, towers, flow dividers in polymer processing 
applications etc. A wealth of information is available on their drag and heat 
transfer characteristics when the liquid phase exhibits the standard Newtonian 
fluid behaviour. In contrast, virtually no analogous results are available even for 
the simplest type of Non-Newtonian fluid behaviour, namely, shear- thinning and 
shear thickening fluids modeled by the usual two-parameter power law model. 

The confined flow of power law fluid past a cylinder with square cross- 
section mounted in a rectangular channel has been investigated. The present 
study is concerned with steady flow at low Reynolds numbers. In order to restrict 
the problem to steady flow, the largest Reynolds number chosen was Re= 50.0 
based on maximum velocity at inlet and the chord length of the square cylinder. 
Two-dimensional computational study has been reported. In the present study for 
the Non-Newtonian fluid, power law model has been used. The power law fluid is 
characterized by two parameters, flow consistency index and power law index. 

The Stream function-Vorticity approach has been employed to solve the 
unsteady Navier-Stocks equations through false transient method. Velocity 
profile and different parameters such as drag coefficient, recirculation length 
have been computed. The effect of power law index (shear thickening 
effect,1.0<« <1.6) and blockage ratio (B = 1/8, 1/6, 1/4) has been investigated. 
Grid independence study also has been studied. 
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Chapter 1 

INTRODUCTION 

1.1 Introduction: 

The incompressible and steady flow of fluids around bluff bodies 
represents an idealization of many industrially important processes, e.g. 
chimneys, tubular heat exchangers, offshore structures and pipelines, 
suspension bridges, towers, masts, wires and flow dividers in polymer 
processing applications etc. The flow past bluff bodies, especially cylinders, has 
been an attraction in all kinds of fluid mechanical investigations for a long time. 
Most of these studies were concerned with the circular cylinder case .In recent 
years; the attention to flow around rectangular cylinders has accelerated, albeit 
only for Newtonian liquids. A wealth of information is available on drag and heat 
transfer characteristics when the liquid phase exhibits the standard Newtonian 
behaviour. In contrast, virtually no corresponding results are available even for 
the simplest type of Non-Newtonian behaviour, namely, shear-thinning and 
shear-thickening modeled by the usual two parameter, power law model fluid. 

1.2 Literature review: 

The literature pertaining to flow past a square cylinder is reviewed below. 
All literatures are concerned with the behaviour of Newtonian fluid. A review of 
the earlier investigations related to the different aspects of the flow over low 
Reynolds numbers has been presented. Viscous flow of Newtonian fluid has 
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been studied extensively in the past. Thus for instance, at low to moderate 
Reynolds numbers, such flows have been computed by Davis and Moore 
(1982), Okajima (1982), Davis et al. (1984), Franke et al. (1990), Okajima et al. 
(1990), Kelkar and Patankar (1992), Mukhopadhyay et al. (1992), Zaki et al. 
(1994), Sohankar et al. (1995), Valencia (1995), Saha et al. (2000), Breuer et 
al. (2000) etc. 

Flow of Newtonian fluid past a square cylinder shows the following 
behaviour with increasing Reynolds number. In the Reynolds number range of 
2-8, the flow is unseparated and symmetric above & below the cylinder with 
respect to the axis aligned with the flow. Once the Reynolds number exceeds 
certain limit, separation at the leading edge is observed. This is followed by 
reattachment over the sides of the cylinder and separation once again at the 
trailing edge and a closed steady recirculation region consisting of two 
symmetric vortices is observed behind the cylinder. The size of recirculation 
region increases with increase in Reynolds number. When a certain threshold 
limit is crossed the vortex shedding can be detected in the wake showing the 
unsteadiness in the flow behind the cylinder, this Reynolds number is called 
critical Reynolds number. 

Many researchers have found different value of critical Reynolds number 
experimentally and numerically. Based on experiment investigations, Okajima 
(1982) found the periodic vortex motion to occur at Reynolds number « 70 thus 
leading to an upper limit of critical Reynolds number <70. A smaller value 
(Re =54) was determined by Kelkar and Patankar (1992) based on the stability 
analysis of the flow. According the numerical study done by Mukhopadhyaya 
(1992), the wake losses its original symmetry at Re =85 and flow becomes 
periodic in the wake region for blockage ratio of 1/4. According to Sohankar et 
al. (1995), the value of critical Reynolds number was found between 50 and 55. 
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Breuer M. et al. (2000) found the unsteadiness at Re =60.0 for the square 
cylinder centrally placed in plane channel for blockage ratio of 1/8. 

As concluded by Okajima (1982) and Franke (1991) in the Reynolds 
number range 100-150, the flow past square cylinders can be considered as 
two-dimensional. In contrast to the circular cylinder flow for which Williamson 
(1996) provides a Reynolds number limit («180) for the onset of three- 
dimensional structures in the wake, no such clear statement can be found for 
the square counterpart. A rough hint has been given by Franke (1991) with 
Re < 300. Therefore, beyond this limit three-dimensional structures are 
expected and subsequently transition to turbulence takes place. So for the 
present study, the steady flow of power law fluids past square cylinder placed in 
rectangular channel, the Re =50.0 was chosen as the upper limit. 

To analyze the flow researchers have used different approaches. For the 
solution of poisson equations for pressure, one compatibility condition has to be 
satisfied. For the staggered-grid this condition is automatically satisfied, 
however this is not the case for the non-staggered grid. Abdallah S. (1987) 
proposed the consistent finite difference method by ensuring the satisfaction of 
compatibility condition on the non-staggered grid. 

The no-slip boundary condition applied to the steady, viscous, 
incompressible flow over a rectangular re-entrant corner boundary leads to 
singularity in the vorticity and pressure fields. Different methods to handle this 
problem proposed by different researchers have been summarized and 
compared by Holstein H. at el. (1981). He shows that discontinuous method 
and Kawaguti method give the good results as compared to others. 
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It is thus safe to conclude that there are no prior results available on the 
confined flow of non-newtonian fluids past a cylinder of square cross-section 
placed in rectangular channel. This work aims to fill this gap in the literature. 


1 .3 Objectives of the present work: 

In the present work, incompressible and steady flow of Non-Newtonian 
fluids (power law model) over a square cylinder mounted in a rectangular 
channel is studied. A finite difference based numerical method Stream function- 
Vorticity scheme has been used to solve the unsteady Navier-Stocks equations 
through false transient methodology. The resulting velocity and pressure fields, 
which are then used to calculate the total, and individual drag coefficients as 
functions of the Reynolds number, power law index and blockage ratio. 
Numerical results for different flow parameters, drag coefficients are presented 
for 5.0<Re<50.0 and 1.0 <n <1.6. The present work is limited for study the 
effect of shear thickening behaviour. The effect of different blockage ratio 
(5 = 1/8, 1/6, 1/4) is also reported. Reynolds number and power law index is 
studied, over a wide range of conditions. A comparison is made with results 
obtained for the Newtonian fluids to that of the Non-Newtonian fluids to 
elucidate the effect of Non-Newtonian fluid characteristics. 
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Chapter 2 

PROBLEM FORMULATION 


This chapter presents the physical problem studied in this work, followed 
by its idealization and finally, the mathematical equations describing the flow 
together with the relevant boundary conditions are presented. The equations are 
non-dimensionalised in order to reduce the number of variables. Also the 
auxiliary equations to calculate the viscous drag coefficient, form drag coefficient, 
non-newtonian viscosity and velocity are presented in this chapter. 


2.1 Problem statement: 

The present work is concerned with the analysis of the steady, laminar 
flow of Non-Newtonian fluid (power law model) past square cylinder mounted in a 
rectangular channel. As mentioned in chapter-1, at low Reynolds number 
(Re<50.0) the flow is two-dimensional, axisymmetric. Figure (2.1) shows the 
schematic diagram of the flow model for the present work. Following dimensions 
(non-dimensionalised with obstacle width) are used: 


LX = 50.0 
17 = 8.0 
4 = 12.5 
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Figure 2.1 : Schematic diagram of flow model 


2.2 Governing equations: 


The governing equations are presented in this section. A rectangular 
coordinate system is chosen for convenience as shown is figure (2.1). The 
velocity components in x & y direction are U , V respectively. The flow is two- 
dimensional and steady. The flow variables do not depend upon the third 
direction. 

Therefore, the velocity components in x-, y- & z- direction are given as 
U'=U'{X',Y') , V'=V\X',Y') and W'=0.0 


Linder these conditions the continuity equation for an incompressible reduces to 


dU' dV 
dX ' + dr 


(2.1) 


And the corresponding x- and y- components of Cauchy equations in stress 

tensor form are written as: 

x-component: 


SIT 

df 


+ U' 


dir ,„eu 


dX' 


+ V- 


ar 


dP' f dr' xx dr' 


ar 


■+ 


+- 


XT 


dX' dY' 


( 2 . 2 ) 
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(2.3) 


y-component: 


fdV 

dV 


dF 


( St 'xy 

dt'yy'] 


+ U ' — + v 

= 

z 


+ — 

v d? 

dX' 

dY'J 

ar 


{ dX' 

dT J 


From the Newton’s law of viscosity for an incompressible fluid the stress can be 
expressed as: 

t - juA (2.4) 


Here, 

ju Coefficient of viscosity 
A=Rate of deformation tensor 


And it is formulated by 

A + 

iJ [a*, sr, J 


(2.5) 


The coefficient of viscosity depends on the local pressure and temperature but 
noton r or A. For the Non-Newtonian fluids the relationship between r&A is 


not the simple proportionality given by equation (2.4) 


For the certain types of the Non-Newtonian fluid, 

t — 7) A (2.6) 

Here, 

tj = Non-Newtonian viscosity 

The assumption of various empirical functions to describe the dependence of 77 
on A (or 7 on r ) corresponds to the various models. 

The simplest type of the model is power law model, in which the non-newtonian 
fluid is characterized by two parameters n,m power law index, flow consistency 
index respectively. 
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Power law model: 


r = < 


m 


|(A:A) 


(«-!)' 


A 


Here, 


a ^=2, t. a , a ., 


If we assume that 


t] = mA 
where, 



Then stress can be represented as: 
r = mAA 


( 2 . 7 ) 


( 2 . 8 ) 

( 2 . 9 ) 


( 2 . 10 ) 


Non-dimensionalisation of equations: 

For the non-dimensionalisation of the velocity terms the maximum velocity 
at inflow plane {U MAX ) is used, for the length terms the obstacle width ( B ) is 

used, for the pressure term pU 2 max is used. 

Following scaling is used for the non-dimensionalisation of equations. 


B B U 


, V = - 


V' 


MAX 


U 


MAX 


w = 


¥' 


U max B 


co = 


co' B 


U 


P = 


P' 


MAX 


p U 2 MAX 


T = ■ 


p U MAX 


A = A' 


f B Y -> 


\Pmax J 
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The definitions of stress tensor from the power law model equation (2.10) 
is substituted in equations (2.3) & (2.4) and after the non-dimensionalisation of 
the resultant equations we get the following equations for the power law fluid. 


x-component: 


dU TT dU TZ dU 3P 1 

+ U + V = + — 

dt 8X dY dX Re 


f d 2 U d 2 U ^ 


+ ■ 


dX 2 dY 2 


+ 2 


dUdA 
dX dX 


■ + 


dU dV_ 
dY + dX 


8A 


dY 


( 2 . 11 ) 


y-component: 

dV TT dV r ,dV dP 1 

H U t- V — = 1 

dt dX dY dY Re 


A 


d 2 v d 2 V 
+ • 


dX 2 dY 2 


_ dV dA f 

+ 2 h 

dY dY 


dU_ dV_)dA 
v dY + dX j 


dX 


( 2 . 12 ) 


Here Reynolds number is defined in following way, 
R 


m 


(2.13) 


By cross-differentiation of the above two equations (2.11) & (2.12) and 
subtraction gives the Vorticity transport equation. 

Vorticity transport equation: 


dco d / \ d ( T \ 1 

— + — (Ua>)+ — (Va>) = — 
dt dX K dY y ’ Re 


dX 2 

4 d 2 A d 2 y/ 


.do) .d 2 co 
A — - + A — r + 


dY 2 


d 2 A d\A 
dX 2 dY 2 


A 


d 2 y/ d 2 y/ 
~dY 2 ~~dX 2 


Re dXdY dXdY 

From the definition of vorticity, stream function elliptic equation is formulated. 
Stream function elliptic equation: 


.dA dco . 

+ 2 + 2 

dX dY 

(2.14) 


d 2 w d 2 w 
— £- + — 4- = -co 
dX 2 dY 2 


(2.15) 

When the equations (2.1 1) & (2.12) are differentiated with respect to x and y 


dA dco 
~dY~dX 
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respectively and added, we get the pressure poisson equation. 


d 2 p d 2 P 
dX 2 + dY 2 


= 4-(Va>)-—(Ua>)+— 
dX 8Y ' Re 


dco dA 
~dX~dY 


dadA dU d 2 A 
dY dX + dX dX 2 


dV d\A 
dY dY 2 


d 2 A 


dU dV 


+ dXdY[dY + dXJ 

(2.16) 


Equations (2.14), (2.15) & (2.16) are solved numerically to get the solution. 


2.3 Initial and boundary conditions: 

The numerical calculations are started with initial guess for all the flow 
field variables over the entire physical domain. Firstly the results were calculated 
for Newtonian fluids (n = 1.0) with different initial guesses and all led to same flow 
field. The best initial guess is that with which the code converges fastest and 
gives the solution. For the Newtonian fluids, fully developed velocity field solution 
gives the fastest convergence. For the Non-Newtonian fluids the calculations are 
started with solution of Newtonian fluids as an initial guess. 

The present problem is solved with stream function-vorticity method. So 
first the boundary conditions are modeled for velocity, and then those boundary 
conditions are formulated for stream function and vorticity. The following 
boundary conditions are used: 

For the momentum equations: 


Inflow plane: Parabolic profile has been used for the velocity. 
U = U(Y)= | , V = 0.0 


' r i") 

y W /2j 
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for the stream function: 



Physically no boundary condition is available for vorticity. So vorticity at 
inflow plane is calculated numerically by equation (2.15). 


• Out flow: There is no unique prescription for the out flow condition. The 
idea is to have specification that would not affect the flow in the 
upstream. So fully developed boundary condition is used. 

V = 0.0 , — = 0.0 
dX 


for the stream function and vorticity: 


dy/ 

~dX 


= 0.0 


dco 

~dX 


= 0.0 


• Confining boundary (channel wall): No-slip boundary condition is used 

U = 0.0 , V = 0.0 

Stream function value at the lower wall: 

n + 1 LY 

w = 

2« + l 2.0 

Stream function value at the upper wall: 

Tf + 1 LY 

y/ = 

2n + l 2.0 

• Obstacle: No- slip boundary condition is used. 

U = 0.0 , V = 0.0 
y/ = 0.0 

No physical boundary condition is available for vorticity at the walls. So the 
vorticity is calculated numerically by equation (2.15) at all the solid walls. 
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For the pressure solution: The pressure boundary conditions are found out 
from the momentum equations. 

• At x=constant boundary 

dP _ Ada 1 f.dU dA dU dA dV dA) 
dX Re dY Rel, dX dX dY dY dX dY J 

• At y=constant boundary 

dP TT Adco 1 f^dVdA dU dA dV dA} 
dY RedX Re l, dY dY dY dX dXdXJ 

2.4 Calculation of Drag: 

After solving equations (2.14), (2.15) and (2.16) numerically to obtain the 
velocity field and pseudo-pressure field, the results are used to estimate the 
individual components of drag coefficient. The drag coefficient represents the 
rate of energy dissipation due to the tangential and normal forces acting on the 
surface of the solid cylinder. Therefore, the total drag comprises of individual 
components, i.e., friction drag due to shear stress and pressure or form drag due 
to pressure. 

The drag coefficient is defined as: 

C D =- — ^ (2.17) 

— pU 2 MAX L 

2 

Here F D is the drag force exerted by the fluid on the cylinder. L represents the 
length in the third direction, as we are solving problem as a two-dimensional, so 
unit length is assumed. Both parts of drag are calculated separately and then 
added. 

Viscous drag coefficient: This part is calculated by integrating the shear force at 
the upper and lower surface of the cylinder 
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(2.18) 



ls: Lower surface 
us: Upper surface 

Form drag coefficient: This part is calculated by integrating the pressure at the 
front and rear side of cylinder. 


\P(Y)dy-\ P(Y)dy 


(2.19) 


1 : Front side of the cylinder 
2: Rear side of the cylinder 

Here P(Y) is the non-dimensional pressure distribution over the cylinder. 


2.5 Calculation of Non-Newtonian viscosity: 


As we have discussed that non-newtonian viscosity can be represented as 
7 = mA 


After the non-dimensionalisation of the equations only 'A' term has to be 
calculated, this is calculated by the following equation: 


A = 


d 2 y/ 

dXdY 


+ 


d 2 y/ d 2 y/ 
~dY 2 ~~dX 2 


/l— 1 


2.6 Calculation of velocity: 


After solving the equations in terms of vorticity and stream function, the 
velocity component are calculated by the following equations: 

dy/ 

~dY 

dy/ 


U = - 


V = — 


dX 
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Chapter 3 


Description of numerical technique 


Fluid flow can be modeled using the principles of conservation of mass 
and Newton’s second law of motion. These laws can be developed in form of 
partial differential equations along with the suitable boundary conditions. These 
differential equations are non-linear and do not admit analytical solutions. 
Numerical approach is the practical alternative for analyzing the flow problem. 

Extensive study has been carried out in the recent decades for developing 
numerical schemes to solve incompressible Navier-Stocks equations in regular 
and multi-dimensional complex geometries. The main difficulty with the 
incompressible flow simulation is the absence of an obvious equation for 
pressure. Specifically the nature of coupling of velocity and pressure variables is 
implicit in nature. When the flow is treated as incompressible, pressure does not 
have the usual thermodynamic meaning. Here it has a 'relative value, which 
adjusts itself instantaneously in such a way that the condition of zero divergence 
is satisfied at all computational grid points. This behaviour is related to the fact 
that the speed of sound becomes infinite in an incompressible fluid. As a 
consequence, the pressure field cannot be calculated by a time-advancement 
procedure. Instead it requires at least one partially implicit determination, which is 
able to take in to account the coupling between the pressure and velocity fields, 
subject to the boundary conditions. 
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3.1 Resolution of pressure and velocity coupling: 

The difficulties associated with the determination of pressure have led to 
the methods that eliminate pressure from the governing equations. Some 
methods are available to solve the incompressible equations. One of the popular 
methods for solving the incompressible Navier-Stocks equations involves 
replacing the primitive variables with stream function and vorticity. If the problem 
is limited to two-dimensional, this method can be applied successfully and easily. 
The elimination of pressure by cross differentiation of two momentum equations 
leads to the vorticity transport equation. This equation when combined with the 
definition of vorticity leads to the basis of stream function-vorticity method. In 
order to determine the pressure at each grid points, it is necessary to solve an 
additional equation, which is referred to as the poisson equation for the pressure. 
Differentiation of equation (2.11) with respect to x and differentiation of equation 
(2.12) with respect to y and add of the results gives the poisson equation for 
pressure. 


3.2 Method of solution: 

The governing equations in each case have been solved on a non- 
staggered grid by using stream function-vorticity method. Since the present work 
is limited to low Reynolds number (5.0<Re<50.0), in this regime the fluid does 
not show any three-dimensionality in the flow, so the present work is limited by 
two-dimensional flow simulation. Among the available methods stream 
function- vorticity method was chosen for the flow simulation because this 
method can be applied easily and successfully for two-dimensional simulation. 

Three equations vorticity transport equation, stream function elliptic 
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equation and pressure poisson equation have to be solved. Since the problem 
is steady, so first the vorticity transport equation and the stream function equation 
are solved simultaneously to get steady state solution for velocity through false 
transient method. Then the steady state velocity solution is used to solve the 
pressure poisson equation for the pressure. The velocity field and pressure field 
are used to find out individual drag coefficients. The flow chart for the complete 
solution is given in figure (3.3) 

3.2.1 Velocity solution: 

The pseudo-unsteady method is used to solve the vorticity transport 
equation with the elliptic stream function equation. In this method the steady 
problem is treated as unsteady problem, the unsteady incompressible equations 
are solved numerically to achieve the steady solutions through convergence. For 
the convective term of the vorticity transport equation, second order upwind 
scheme has been used. Diffusion terms of the vorticity transport equation have 
been approximated by second order central difference. The same second order 
difference scheme have been employed for the second order derivative in the 
stream function elliptic equation 

Following algorithm is used to get velocity solution: 

• Specify initial values for stream function and vorticity at time t = 0.0 

• Solve the vorticity transport equation for vorticity at each interior grid 
points at t + At 

• Iterate the new stream function value at all points by solving the elliptic 
equation using new vorticity values at interior points. 

• Find the velocity components from 

dY dX 
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• Determine values of vorticity on the boundaries using new value of 

stream function and vorticity at interior points. 

• Return to step 2 if the solution is not converged 

3.2.2 Pressure solution: 

Solution for the pressure is obtained from the poisson equation 
(divergence of the momentum equation) with the Neumann boundary conditions. 
Solution for the poisson equation exists only if a compatibility condition is 
satisfied. This condition relates the source of poisson equation (2.16) and the 

Neumann boundary conditions (Green’s theorem). 

1 1 

J J cr dxdy = \P nb dS 

x=0 y = 0 

Here, 

cr : Source term of the poisson equation (2.16) 

«6:Outward normal to boundary contour S, enclosing the solution 
domain 

This compatibility condition is not automatically satisfied on the non- 
staggered grid. Failure to satisfy the compatibility condition leads to non- 
convergent iterative solutions. A number of methods are available by ensuring 
the convergence. In the present work a method developed by Abdallah S. (1986) 
has been used. Consistent finite difference approximations for the pressure 
equation with Neumann boundary conditions have been developed to exactly 
satisfy the compatibility. He developed the scheme for the equations for the 
Newtonian fluids. That scheme has been modified for the Non-Newtonian fluids 
(power law model). Here some outlines of that method for the Newtonian fluids 
with the specified grid arrangement are presented. 
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Consistent Finite Difference approximations: Following finite difference grid 
arrangement is used for this scheme. 
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Figure 3.1: Grid arrangement for solution of pressure poisson equation 

Pressure poision equation for the Newtonian fluid is written as: 

Pxx + Pyy = — 

The subscripts x, y refer to the partial derivatives with respect to x, y. 

The equation is approximated at grid points (i, j) using second order accurate 
formulas in the following way: 

to - ra\ - to - ra). + to + Ua). - (. P r + Ua>\ =0.0 
e,w,n,s are the mid points on the lines joining the grid points. 

The pressure gradients P x and P y are approximated using second order 
accurate formulas at the locations e,w,n,s as 
(Px)„ =(?,.! -P, 

to),=to- j v.) /A 

to).=to..-' p u)^ 

h : grid spacing (uniform in both direction) 
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The values of velocity and vorticity at e,w,n,s are calculated from their values at 
the grid locations by averaging 

For consistency, the pressure gradients P x is evaluated at the location x = h / 2 
and not at x = 0.0. All the pressure gradients are evaluated not at the exact 
boundary grid points, but at the point that is middle of boundary grid point and the 
next grid point to boundary grid point in the domain. By using this approximation, 
the compatibility condition is satisfied. 


3.3 Grid system: 

The computational domain is divided into set of rectangular cells and a 
non-staggered arrangement is used such that all variables are defined at the 
corners of the cells. Grid points exactly lies on the cylinder surface. Two- 
dimensional simulation has been carried out on uniform grids and grid 
independence study has been discussed in chapter-4. 


3.4 Finite Difference Discretization: 

The present problem has been solved with the finite difference based 
methodology. In this method, the derivatives in the partial differential equations 
are replaced by finite difference quotients derived from the Taylor series 
expression. The time derivative is discritized explicitly. For the convective term 
of the vorticity transport equation, second order upwind scheme has been used. 
Diffusion terms of the vorticity transport equation have been approximated by 
second order central difference. 
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Second order upwind scheme: The following Discretization (sample) is used in 
the second order upwind scheme. 


ij+1 
-I- w 


i-U 


i+lj 


4s 
i.j-1 

Figure 3.2: Grid arrangement for the second order upwind scheme 


d A., \ _ U R CO R U 


Here, 






DX 


U L = U ‘ +U -> 

L. r\ 


if U R > 0.0 then co R = co i 
if U R < 0.0 then co R - a> M 


if U L > 0.0 then co L = a>._ x 
if U L < 0.0 then co L = co i 


The discritization of vorticity transport equation (2.14) is presented here. 
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3.5 Treatment of corner points: 

The comer point treatment is very much important. It is well known that 
the no-slip, plane, viscous, incompressible fluid flow over a rectangular re-entrant 
corner boundary leads to a singularity in the vorticity and the pressure fields. The 
presence of the singularity implies that the stream function and vorticity do not 
have Taylor series expansion in the regions that include corners. Consequently 
the usual discritization methods cannot be applied. Two methods Kawaguti 
method and discontinuous method were tried and there were no so much 
effective variations in the results. The results, which are reported here, are based 
on discontinuous method. In the discontinuous method, at the corner point two 
values are calculated for the vorticity and pressure obtained from the both (x- and 
y-) direction equations. 
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3.6 Numerical stability: 


For accuracy, the mesh size must be small enough to resolve the 
expected spatial variation in all the dependent variables. But the usage of very 
fine grid makes the calculation very expensive. So on the basis of grid 
independence test, a proper mesh size is chosen. Once the mesh has been 
chosen, the condition on the time step necessary to prevent numerical 
instabilities are determined from the combination of CFL (Courant-Friedriches- 
Law) condition and the restriction on the grid Fourier number. According to CFL 
condition, the distance the fluid travels in one time increment must be less then 
one space increment. This leads to a constraint on the time step in the form 


At < 


AX AY 1 

m ’ n i 


When the viscous diffusion terms are important, the condition necessary to 
ensure stability is dictated by the restriction on the grid Fourier number and can 
be shown to be 


At < 


RE 

2.0 


AX 2 A y 2 
AX 1 + AY 2 * 


The final time step is the minimum of the values obtained from. 
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Figure 3.3: Flow chart for the solution of problem 
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Chapter-4 


Results and Discussion 


The flow of Non-Newtonian fluids (power law model) past a square 
cylinder placed in a rectangular channel has been analyzed in the present work. 
The present study has been done in the low Reynolds number range 
(5.0 < Re <50.0). The streamlines patterns, iso-vorticity patterns, recirculation 
length, individual coefficient of drags have been computed numerically. The 
present study is concerned with the shear-thickening fluids or dilantant fluids 
(n >1.0). In the present work, study has been done for 1.0 <n <1.6. The effect of 
the blockage ratio (B = 1/8, 1/6, 1/4) has also been reported. The numerical 
results are presented and discussed in this chapter. 

Time-step used in the problem first was calculated as described in section 
3.6.The final time-step used in the problem varies with the different parameters 
such as Reynolds number, power law index, and blockage ratio. For 5 = 1/8, 
n =1.0 and Re =50.0 case, Ar(time step)=0.001 was used. For the lower value 
of Reynolds number (Re < 50.0) and for n > 1 .0, further lower value of the time- 
step was used as compared to upper given value. 


4.1 Grid Independence Test: 

Three different spatial grid spacing (uniform and same in both direction) 
(Ax = 0.2,0.1,0.05) were used and results were compared for the one case. The 
following figure 4.1 shows the effect of the different spatial grid spacing. 
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Figure 4.1 Grid independence test (Re = 20.0,5 = l/8,« = 1.0) 


There is not much difference between results obtained with 
Ax = Ay = 0.1(501x81) & Ax = Ay = 0.05(1001x161), but the time taken by the code 
with finest grid Ax = Ay = 0.05 (1001x1 61) was very much as compared to that of 
Ax = Ay = 0.1(501x81), so finally Ax = Ay = 0.1 (501x81) uniform grad spacing 
was chosen for solving the present problem numerically. 


4.2 Streamlines and Iso-vorticity Plots: 

Streamlines and iso-vorticity plots provide useful information about the 
detailed kinematics of the fluid flow. These plots provide the qualitative 
information about the flow. Streamlines contour plots show the separation of the 
flow, reattachment in the flow and recirculation region. Some representive stream 
function and vorticity plots have been presented herein. Different plots have been 
categorized in the following way: 
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4.2.1 Effect of the blockage ratio: 

Figure 4.2 to 4.4 shows the streamlines patterns corresponding to 

5 = 1/8, 1/6, 1/4 at » = 1.0 and Re = 50.0. 

4.2.2 Effect of the Reynolds number: 

In this section streamlines patterns figure 4.5 to 4.7 have been drawn to 
show the effect of Reynolds number at constant value of power law index n = 1.2 

6 blockage ratio 5 = 1/8. It is clear from figure 4.4 to 4.7 (corresponding to 

Re = 10.0,20-0,50.0 Respectively) that as the Reynolds number increases; the 
length of recirculation region increases. 

4.2.3 Effect of the power law index: 

Figure 4.8 to 4.11 shows the streamlines patterns corresponding to 
n = 1.0,1 .2,1 .4,1 .6 at Re = 50.0 and 5 = 1/6. These plots have been drawn to show 
the effect of power law index. Similarly figure 4.12 to 4.15 shows the iso-vorticity 
plots corresponding to n = 1.0,1 .2,1 .4,1 .6 respectively for 5 = 1/6 and Re = 50.0. 


4.3 Recirculation Length: 

For the flow problems at steady state, the recirculation length is also 
important parameter. Recirculation length for the square cylinder has been 
computed numerically. This parameter has been computed as a function of 
Reynolds number, power law index, and blockage ratio. Some representive 
values have been summarized in the following tables to show the effect of all the 
three parameters (Re,n,5) 
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4.3.1 Effect of Reynolds number: 

The following table 4.1 shows the recirculation length as a function of 
Reynolds number at constant value of blockage ratio and power law index. 

Table 4.1 : Effect of the Reynolds number on the recirculation length 


Reynolds Number 

Recirculation Length 

(5 = 1/8&b = 1.2) 

5.0 

0.23 

10.0 

0.50 

20.0 

1.00 

50.0 

2.72 


It’s clear for the above table that as the Reynolds number increases, the 
recirculation length increases. 

4.3.2 Effect of Power law index: 

The following table 4.2 shows the recirculation length as a function of 
power law index at constant value of Reynolds number and blockage ratio. 

Table 4.2: Effect of power law index on the recirculation length 


Power law index 

Recirculation length 

(Re = 50.0&5 = 1/8) 

1.0 

2.51 

1.2 

2.72 

1.4 

2.91 

1.6 

3.05 


It’s clear from the table 4.2 that as the value of power law index increases, 
the recirculation length also increases. But this happens only at lower blockage 
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ratio (5 = 1/ 8,1/6). It was found that at B = 1/4, the recirculation length does not 
very with respect to power law index much. 


4.3.3 Effect of the blockage ratio: 

The following table 4.3 shows the effect of blockage ratio on the 
recirculation length at constant value of Reynolds number and power law index. 
Table 4.3: Effect of blockage ratio on recirculation length 


Blockage Ratio 

Recirculation Length 

(Re = 50.0 &n =1.2) 

1/8 

2.72 

i/6 

2.01 

1/4 

1.20 


It’s clear for the table 4.3 that recirculation length decreases with increase in 
blockage ratio. 


4.4 Coefficient of drag: 

One of the important characteristic quantities of the flow around the 
cylinder is the drag coefficient. While the streamlines and iso-vorticity patterns 
provide information about the kinematics of flow, the drag phenomena provide 
the overall hydrodynamic picture of the flow. In the region of small Reynolds 
number the coefficient of drag varies strongly with Reynolds number. The 
contribution of the viscous and pressure forces to the total force is of the same 
magnitude. The coefficient of drag has been computed with various parameters. 
The values of the individual coefficient of drags and the total coefficient of drag 
are summarized in the following tables with respect to different parameters. 
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Table 4.4: Viscous, form and total drag coefficient (B = 1/8 , n = 1.0) 


Reynolds 

Number 

Viscous Drag 
Coefficient ( C D ) 

Form Drag 
Coefficient ( C D ) 

Total Drag 

Coefficient ( C D ) 

5.0 

1.860 

5.225 

7.085 

10.0 

1.019 

3.272 

4.291 

20.0 

0.549 

2.249 

2.798 

50.0 

0.219 

1.588 

1.807 


Table 4.5: Viscous, form and total drag coefficient (B = l/8,n = 1.2 ) 


Reynolds 

Number 

Viscous Drag 

Coefficient ( C D ) 

Form Drag 

Coefficient ( C D ) 

Total Drag 

Coefficient ( C D ) 

5.0 

1.699 

2.659 

4.348 

10.0 

0.958 

1.823 

2.781 

20.0 

0.536 

1.411 

1.947 

50.0 

0.227 

1.210 

1.427 


Table 4.6: Viscous, form and total drag coefficient (5 = 1/8,« = 1.4) 


Reynolds 

Number 

Viscous Drag 

Coefficient ( C D ) 

Form Drag 

Coefficient ( C D ) 

Total Drag 

Coefficient ( C D ) 

5.0 

1.618 

1.838 

3.456 

10.0 

0.944 

1.299 

2.243 

20.0 

0.547 

1.058 

1.605 

50.0 

0.245 

0.99 

1.235 


. 29 




















































Table 4.7: Viscous, form and total drag coefficient (i? = l/8,« = 1.6) 


Reynolds 

Number 

Viscous Drag 

Coefficient ( C D ) 

Form Drag 

Coefficient ( C D ) 

Total Drag 

Coefficient (C D ) 

5.0 

1.518 

0.962 

2.480 

10.0 1 

0.924 

0.737 

1.661 

20.0 

0.555 

0.666 

1.221 

50.0 

0.264 

0.704 

0.968 


Table 4.8: Viscous, form and total drag coefficient (B = l/6,n=l.O) 


Reynolds 

Number 

Viscous Drag 
Coefficient ( C Df ) 

Form Drag 

Coefficient ( C D ) 

Total Drag 

Coefficient ( C D ) 

5.0 

2.242 

6.415 

8.657 

10.0 


3.825 

5.004 

20.0 

0.618 

2.533 

3.151 

50.0 

0.244 

1.743 



Table 4.9: Viscous, form and total drag coefficient (B = l/6,n = 1.2) 


Reynolds 

Number 

Viscous Drag 

Coefficient ( C D ) 

Form Drag 

Coefficient ( C Dp ) 

Total Drag 

Coefficient ( C D ) 

5.0 

2.162 

3.386 

5.448 

10.0 

1 

1.155 

2.149 

3.304 

20.0 

0.621 

1.595 

2.216 

50.0 

0.260 

1.315 

1.575 
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Table 4.10: Viscous, form and total drag coefficient (B = 1/6 , n = 1.4) 


Reynolds 

Number 

Viscous Drag 

Coefficient ( C D ) 

Form Drag 

Coefficient ( C D ) 

Total Drag 

Coefficient ( C D ) 

5.0 

2.175 

2.353 

4.528 

10.0 

1.184 

1.536 

2.72 

20.0 

0.655 

1.188 

1.843 

50.0 

0.286 

1.070 

1.356 


T able 4.11: Viscous, form and total drag coefficient ( B = 1 / 6, n = 1 .6 ) 


Reynolds 

Number 

Viscous Drag 

Coefficient ( C D ) 

Form Drag 

Coefficient ( C D ) 

Total Drag 

Coefficient ( C D ) 

5.0 

2.169 

1.262 

3.421 

10.0 

1.202 

0.862 

2.064 

20.0 

0.680 

0.817 

1.497 

50.0 

0.309 

0.754 

1.063 


Table 4.12: Viscous, form and total drag coefficient (B = l/A,n = 1.0) 


Reynolds 

Number 

Viscous Drag 

Coefficient ( C D ) 

Form Drag 

Coefficient ( C D ) 

Total Drag 

Coefficient ( C D ) 

5.0 

3.396 

10.190 

13.586 

10.0 

1.702 

5.624 

7.326 

20.0 

0.854 

3.437 

4.291 

50.0 

0.332 

2.190 

2.512 
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Table 4.13: Viscous, form and total drag coefficient (B = 1/ 4,n = \.2) 


Reynolds 

Number 

Viscous Drag 
Coefficient ( C D ) 

Form Drag 

Coefficient ( C D ) 

Total Drag 

Coefficient ( C D ) 

5.0 

3.695 

5.956 

9.651 

10.0 

1.858 

3.415 

5.273 

20.0 

0.943 

2.323 

3.175 

50.0 

0.380 

1.643 

2.123 


Table 4.14: Viscous, form and total drag coefficient (5 = l/4,n = 1.4) 


Reynolds 

Number 

Viscous Drag 

Coefficient ( C D ) 

Form Drag 

Coefficient (C Dp ) 

Total Drag 

Coefficient ( C D ) 

5.0 

4.170 

4.503 

8.673 

10.0 

2.109 

2.602 

4.711 

20.0 

1.083 

1.704 

2.787 

50.0 

0.447 

1.324 

1.771 


Table 4.15: Viscous, form and total drag coefficient (B = l/4,« = 1.6) 


Reynolds 

Number 

Viscous Drag 

Coefficient ( C D ) 

Form Drag 

Coefficient ( C Dp ) 

Total Drag 

Coefficient ( C D ) 

5.0 

4.667 

2.661 

7.328 

10.0 

2.369 

1.553 

3.992 

20.0 

1.226 

1.023 

2.249 

50.0 

0.516 

0.890 

1.406 
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The individual coefficient of drags and total coefficient of drag have been 
presented in the tabular form for the various parameters. In the next session 
coefficient of drag will be represented in the graphical form to show the effect of 
Reynolds number, blockage ratio and power law index. 


4.4.1 Effect of Reynolds number and Power law index: 

Figure 4.16 to 4.18 have been plotted for Coefficient of Drag v/s Reynolds 
number for three different blockage ratio (1/8, 1/6, 1/4). In all the three graphs, the 
plots have been shown for different values of the power law index 
(1 .0,1 .2,1 .4,1 .6). It’s clear from the graph that as the Reynolds number increases 
by keeping blockage ratio and power law index constant, the coefficient of drag 
decreases. Similarly at the constant value of blockage ratio and Reynolds 
number, the coefficient of drag decreases with increase in the power law index or 
increasing the shear thickening effect. 

4.4.2 Effect of Reynolds number and blockage ratio: 

Figure 4.19 to 4. 22 have been plotted for coefficient of drag v/s Reynolds 
number for different values of power law index (n = 1.0, 1.2, 1.4, 1.6). In all 4 graphs, 
the effect of blockage ratio is also presented. It’s clear form the graphs that as 
the blockage ratio increases, the coefficient of drag also increases. 


33 



Figure 4.2 


Streamlines pattern 

Blockage Ratio: 1/8, Reynolds Number: 50.0 
Power Law Index: 1.0 
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Figure 4.3 


Streamlines pattern 
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Figure 4.4 


Streamlines pattern 

Blockage Ratio: 1/4, Reynolds Number. 50.0 
Power Law Index: 1.0 
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Figure 4.5 


Streamlines pattern 

Reynolds Number: 10.0, Blockage Ratio: 1/8 
Power Law Index: 1.2 


37 




Figure 4.6 


Streamlines pattern 

Reynolds Number 20.0, Blockage Ralio:l/8 

power Law Index. L2 




Figure 4.7 


Streamlines pattern 

Reynolds Number. 50.0, Blockage Ratio: 1/8 
Power Law Index: 1.2 
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Figure 4.8 


^ucdinniies pattern 


dockage R^fo- 1/6 '°’ Reynolds Numb <=n 5C 
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Figure 4.9 


Streamlines pattern 

Power Law Index: 1.2, Reynolds Number: 50.0, 
Blockage Ratio: 1/6 
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Figure 4. 10 

Streamlines pattern 

Po wer Law Index: 1 . 4 , Reynolds Number 50 0 
Blockage Ratio: 1/6 U U ' 


A"> 


Figure 4.11 


Streamlines pattern 

Power Law Index: 1.6 , Reynolds Number: 50.0, 
Blockage Ratio: 1/6 




Iso-Vortieity Plot 
Power Law lndex=1 .0 

Blockage Ratio: 1/6, Reynolds Number: 50.0 


Figure 4.1 3 


Iso-Vorticity Plot 
Power Law lndex=1 .2 

Blockage Ratio: 1/6, Reynolds Number: 50.0 
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Figure 4.1 4 


Iso-Vortieity Plot 
Power Law lndex=1 .4 

Blockage Ratio: 1/6, Reynolds Number: 50.0 
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Figure 4. 15 


Iso-Vorticity Plot 
Power Law lndex=1 .6 

Blockage Ratio: 1/6, Reynolds Number: 50.0 
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Total Drag Coefficient ~ 


Effect of Reynolds Number & Power Law Index 
Bl ockag e Rati o: MB 



Figure 4.16 



Total Drag Coefficient ~ 


Effect of Reynolds Number & Power Law Index 
Blockage Ratio: 1/6 



Figure 4.1 7 




Effect of Reynolds Number & Power Law Index 
Blockage Ratio: 1/4 
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Total Drag Coefficient - 


Effect of Reynolds Number & Blockage Ratio 
Power Law Index : 1 .0 
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Figure 4.20 



Total Drag Coefficient ; 


Effect of Reynolds Number & Blockage Ratio 
Power Law Index : 1 .4 




Total Drag Coefficient - 


Effect of Reynolds Number & Blockage Ratio 
Power Law Index : 1 .6 




Chapter 5 


CONCLUSION AND SUGGESSIONS FOR 
FUTURE WORK 


In the present work, Non-Newtonian (power law model) fluids flow past 
square cylinder placed in rectangular channel has been analyzed. The present 
problem has been solved with Stream function-Vorticity method. Problem has 
been solved on the uniform grid spacing. In particular, the detailed kinematics of 
the flow in terms of the values of stream function and vorticity has been obtained. 
Further, the values of stream function and vorticity were used to calculate 
pseudo-pressure field. Finally velocity field and pressure field were used to 
calculate the individual coefficient of drags and total coefficient of drag. The 
present study covers the low Reynolds number range (5.0 < Re < 50.0) . The 
shear thickening effect (1.0 <n <1.6) and blockage ratio (1 / 8,1 / 6,1 / 4) has been 
studied. The following points can be concluded from the present work 

1. The flow remains steady below 50.0 Reynolds number for all values of 
power law index. In the wake of the cylinder, a steady axi-symmetric 
recirculation region is observed. The length of this recirculation region 
increases with increase in Reynolds number. 

2. At a constant value of power law index and blockage ratio, the viscous 
drag coefficient, form drag coefficient and total drag coefficient decrease 

with increase in Reynolds number. 



3. At a constant valcie of Reynolds number and blockage ratio, the length of 
recirculation region increases with increase in power law index. This can be 
explained by analyzing the figure 5.1. If is clear from the figure that as the 
power law index increases (implying high viscosity near the cylinder surface) 
the stream function contour line of the same value (\j/ =0.1) goes away from 
the cylinder. So this creates long recirculation region in the wake. The total 
coefficient of drag decreases with increase in power law index. 

4. At a constant value of Reynolds number and power law index, the viscous 
drag coefficient, form drag coefficient and total drag coefficient increase with 
increase in blockage ratio. 

Suggestion for future work: 

The work reported herein should be extended to: 

1 .elucidate the role of shear thinning on drag as well as on the streamlines 
and vorficify plots. 

2.calculafe the heat transfer characteristics such as Nusself number for 
shear thickening and shear thinning fluids. 


Stream lines for different value of Power Law Index for 
same value of stream function (^-0.1 ) 

Re=50.0 , B=1/6 

Figure 5.1 
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